201605.00514v1 


chinaXiv 


AL CHANGE BIOLOGY 


BIOENERGY 


GCB Bioenergy (2016), doi: 10.1111/gcbb.12316 


Ecosystem model parameterization and adaptation for 
sustainable cellulosic biofuel landscape design 


JOHN L. FIELD'?, ERNIE MARX!, MARK EASTER', PAULR. ADLER? and 

KEITH PAUSTIAN!4 

1Natural Resource Ecology Laboratory, Colorado State University, 1499 Campus Delivery, Fort Collins, CO 80523, USA, 
Department of Mechanical Engineering, Colorado State University, 1374 Campus Delivery, Fort Collins, CO 80523, USA, 
Pasture Systems and Watershed Management Research Unit, United States Department of Agriculture-Agricultural Research 
Service, University Park, PA 16802, USA, “Department of Soil and Crop Science, Colorado State University, 1170 Campus 
Delivery, Fort Collins, CO 80523, USA 


Abstract 


Renewable fuel standards in the US and elsewhere mandate the production of large quantities of cellulosic 
biofuels with low greenhouse gas (GHG) footprints, a requirement which will likely entail extensive cultiva- 
tion of dedicated bioenergy feedstock crops on marginal agricultural lands. Performance data for such systems 
is sparse, and non-linear interactions between the feedstock species, agronomic management intensity, and 
underlying soil and land characteristics complicate the development of sustainable landscape design strategies 
for low-impact commercial-scale feedstock production. Process-based ecosystem models are valuable for 
extrapolating field trial results and making predictions of productivity and associated environmental impacts 
that integrate the effects of spatially variable environmental factors across diverse production landscapes. 
However, there are few examples of ecosystem model parameterization against field trials on both prime and 
marginal lands or of conducting landscape-scale analyses at sufficient resolution to capture interactions 
between soil type, land use, and management intensity. In this work we used a data-diverse, multi-criteria 
approach to parameterize and validate the DayCent biogeochemistry model for upland and lowland switch- 
grass using data on yields, soil carbon changes, and soil nitrous oxide emissions from US field trials spanning 
a range of climates, soil types, and management conditions. We then conducted a high-resolution case study 
analysis of a real-world cellulosic biofuel landscape in Kansas in order to estimate feedstock production poten- 
tial and associated direct biogenic GHG emissions footprint. Our results suggest that switchgrass yields and 
emissions balance can vary greatly across a landscape large enough to supply a biorefinery in response to 
variations in soil type and land-use history, but that within a given land base both of these performance fac- 
tors can be widely modulated by changing management intensity. This in turn implies a large sustainable cel- 
lulosic biofuel landscape design space within which a system can be optimized to meet economic or 
environmental objectives. 
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environmental impacts and causing minimal disruption 


I ; 
mboduetion to agricultural commodity markets, which could lead to 


The Energy Independence and Security Act of 2007 
expanded the US renewable fuel standard to require the 
use of large quantities of ‘advanced’ and ‘cellulosic’ bio- 
fuels with lifecycle greenhouse gas (GHG) emissions 
reductions of 50% and 60%, respectively, relative to a 
conventional gasoline baseline (110th Congress of the 
United States, 2007). This mandate is predicated on the 
wide availability of biomass feedstocks with low direct 
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indirect leakage effects (Searchinger et al., 2008). A vari- 
ety of low-impact feedstock provisioning strategies have 
been envisioned including the collection of agricultural 
residues, forestry residues, and municipal wastes, as 
well as the cultivation of dedicated woody and herba- 
ceous crops on marginal or non-agricultural lands 
(Campbell et al., 2008; Robertson et al., 2008; Tilman 
et al., 2009). In particular, perennial grasses such as 
switchgrass, Miscanthus, bioenergy-optimized sugar- 
cane, and mixed prairie species have been identified as 
promising cellulosic feedstock crops due to their high 
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yield potential, relatively low input requirements, high 
nitrogen use efficiency, and ability to sequester carbon 
through soil organic matter formation (Vogel etal., 
2002; Tilman et al., 2006; Heaton et al., 2008; Walter 
et al., 2014). A comprehensive review of US feedstock 
availability by the Department of Energy suggests that 
widespread cultivation of high-yielding dedicated 
perennial grasses will be necessary to achieve the target 
of displacing 30% of U.S. petroleum consumption with 
biofuels (U.S. Department of Energy, 2011), and it has 
been estimated that biofuel supply chains based on such 
feedstocks will have highly favorable lifecycle GHG 
impacts as compared to first-generation biofuel tech- 
nologies (Schmer et al., 2008; Davis et al., 2012; Wang 
et al., 2012). 

Despite their promise, agronomic experience with 
perennial grass feedstock crops is still relatively limited, 
and questions around the best management practices to 
balance the often-competing goals of maintaining high 
yields while minimizing environmental impacts have 
not yet been resolved. The debate between the relative 
merits of low-intensity cultivation over large areas 
(often referred to as ‘land-sharing’) vs. more intensive 
production on a more limited land base (‘land-sparing’) 
is still being waged (Anderson-Teixeira et al., 2012). Pro- 
cess-based ecosystem modeling can play an important 
role in extrapolating limited existing perennial grass 
field trial results to make more general estimates of 
productivity, environmental impacts, and optimal man- 
agement and landscape design strategies. 


Modeling cellulosic feedstock yields & environmental 
impacts 


The use of crop models for assessing management- 
environment interactions and predicting bioenergy feed- 
stock crop productivity was thoroughly reviewed by 
Nair et al. (2012). Crop models such as APSIM, BioCro, 
and ALMANAC have been applied at regional or 
national scales to assess the productivity of first- and 
second-generation feedstock crops as affected by broad- 
scale climate-soil associations (Bryan et al., 2010; Miguez 
et al., 2012; Behrman et al., 2014). Bioenergy system 
design is not just a question of yield, however, and 
understanding the biogeochemical cycling of carbon, 
nitrogen, and water through these systems is essential 
for quantitatively evaluating their sustainability (Robert- 
son etal., 2011). Perennial grass feedstock crops are 
often associated with a high potential for soil carbon 
sequestration. One recent meta-analysis suggests that 
switchgrass increases soil organic carbon (SOC) levels at 
a median rate of ca. 0.7 t C ha! yr! when cultivated 
on carbon-depleted agricultural soils, though perfor- 
mance is more neutral on pastureland or areas that 


were not previously cropped (Qin et al., 2015a). Appli- 
cation of nitrogen fertilizers is typically required to 
replace losses during harvest and maintain yield levels, 
but emissions of nitrous oxide (N2O, a byproduct of soil 
microbial activity and a potent GHG) increase non-line- 
arly with increasing N rate (Hoben et al., 2011; Shcher- 
bak et al., 2014) and can threaten the overall lifecycle 
GHG footprint of any bioenergy system based on feed- 
stocks with inefficient nitrogen cycling (Crutzen et al., 
2008). The biogeochemical cycling of C, N, and HO 
and associated fluxes of biogenic GHGs are tightly 
linked in all agroecosystems by fundamental mecha- 
nisms including plant tissue stoichiometry, photosyn- 
thetic pathway (C3; vs. C4 vs. CAM), stomatal 
conductance, and microbial mineralization/immobiliza- 
tion, and are sensitive to local climate, soil type, and 
land-use history. 

Detailed reviews of process-based biogeochemical 
model use to capture these interactions in the context of 
bioenergy system sustainability assessment are pro- 
vided by Thomas ef al. (2013) and Robertson et al. 
(2015). The CENTURY model was among the first to be 
applied to bioenergy sustainability assessment, and it 
and its derivative DayCent model have been widely 
used to evaluate corn grain production, corn stover 
removal, and the dedicated cultivation of switchgrass 
and Miscanthus from the level of individual sites to 
national scales (Sheehan et al., 2003; Kim & Dale, 2005; 
Chamberlain et al., 2011; Davis et al., 2012; Lee et al., 
2012; Duval et al., 2013). The Environmental Policy Inte- 
grated Climate model has also been applied extensively 
to bioenergy feedstocks in the context of economic anal- 
yses (Jain et al., 2010; Egbendewe-Mondzozo et al., 2011) 
and environmental sustainability assessments (Gelfand 
et al., 2013) at scales from regional (Zhang et al., 2010) to 
global (Kang et al., 2014). 

Experience with using models to assess the produc- 
tivity and biogeochemical implications of fine-scale vari- 
ations in environmental factors such as soil type or 
topography across landscapes or individual farms is 
more limited. Studies of switchgrass cultivation in the 
southeastern US using crop production models have 
come to different conclusions as to whether yield is 
sensitive (Woli, 2012) or insensitive (Persson et al., 2011) 
to underlying soil type. At the finest spatial scale, a 
group at Idaho National Laboratory has used corn yield 
data from a precision agriculture system in Iowa to 
drive DayCent and models of soil erosion and identify 
areas of low nitrogen use efficiency to target for conver- 
sion to switchgrass (Abodeely et al., 2013). Many studies 
conduct large numbers of fine-scale simulations to make 
regional-scale estimates of feedstock productivity and 
environmental performance, though they typically do 
not report on soil-climate-management interactions 
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explicitly but rather emphasize more aggregate descrip- 
tions of landscape performance (e.g. Davis et al., 2012; 
Gelfand et al., 2013). 


Challenges in bioenergy ecosystem modeling 


The use of ecosystem models to assess different bioen- 
ergy landscape design strategies is complicated by 
challenges around the representation of marginal 
lands, adequate bioenergy crop parameterization, and 
selection of the most appropriate spatial resolution 
and agronomic management practices for simulation. 
While biogeochemical process models are increasingly 
used to simulate conversion of marginal agricultural 
lands to bioenergy feedstock cultivation (Bandaru 
et al., 2013; Qin et al., 2015b), these scenarios are par- 
ticularly challenging from a modeling perspective. The 
definition of ‘marginal’ land itself is not straightfor- 
ward or consistent across studies, and the term is 
applied to both land with low productivity potential 
and to land with ample productivity but vulnerable to 
long-term degradation (erosion, loss of soil organic 
matter, etc.) under conventional cropping systems. In 
many bioenergy assessment studies the marginal des- 
ignation has been based on unfavorable biophysical 
properties as judged using land suitability ratings 
(Gelfand et al., 2013) or remote sensing techniques (Cai 
et al., 2011). Another more direct basis for the designa- 
tion considers past transitions in an out of agricultural 
production as inferred from land-use datasets (Camp- 
bell et al., 2008), remote sensing (Wright & Wimberly, 
2013), or sector-level economic modeling (Swinton 
et al., 2011). When considering lands designated as 
marginal based on their low productivity potential, 
bioenergy feedstock crops are often not immune to the 
factors that make such lands challenging for conven- 
tional crops. Recent perennial grass field trials pur- 
posely conducted on marginal sites indicate reduced 
productivity relative to performance on the prime 
lands typically encountered at many agricultural field 
stations (Mooney et al., 2009; Shield et al., 2012; Boyer 
et al., 2013). From an ecosystem modeling perspective, 
accurate assessment is only possible to the extent that 
underlying biophysical limitations on productivity (e.g. 
unfavorable climates, soil texture extremes, shallow 
soils, low soil organic matter levels, site drainage 
problems, slope, vulnerability to erosion, etc.) are rep- 
resented directly or indirectly in model data inputs 
and the processes simulated. 

It is well understood that process-based ecosystem 
models require proper parameterization specific to the 
agroecosystems being simulated in order to achieve 
good performance, and that such models have limited 
predictive power when extrapolated far beyond their 


parameterization scope (Thomas et al., 2013). However, 
many bioenergy assessment studies are based on mod- 
els parameterized under prime conditions and then 
extrapolated to highly marginal sites, or lacking an 
explicit independent validation of performance (e.g. 
Gelfand et al., 2013; Kang et al., 2014). In the case of the 
DayCent model, parameterization typically involves 
adjusting study site and crop parameters by hand in 
order to match observed real-world performance for a 
small number of field trial cases for which extensive 
data are available (Del Grosso et al., 2011). While this 
approach can often yield a high degree of fidelity across 
a range of performance criteria for the sites in question 
(Hudiburg et al., 2015), the very large number of empiri- 
cally-determined crop and site parameters in the model 
makes the process vulnerable to over-parameterization. 
In such cases, model fit to the training dataset is 
improved via mechanisms that lack broader underlying 
ecological significance, reducing the generality of the 
resulting model for other geographic areas, environ- 
mental conditions, or management regimes (Necpalova 
et al., 2015). It is also possible to introduce bias with the 
selection of the parameterization cases themselves, if 
the researcher gravitates to focus on studies that 
confirm a priori assumptions of how a system ‘should’ 
perform. 

There are additional challenges around the spatial 
resolution of landscape modeling and assumptions 
about crop agronomic management. Management fac- 
tors including tillage intensity (Adler etal., 2007), 
fertilizer application rate (Davis et al., 2013), and rota- 
tion length (Pyörälä et al., 2014) can potentially change 
the lifecycle GHG performance of a bioenergy system 
from positive to negative, an effect termed the ‘man- 
agement swing potential’ (Davis et al., 2013). Manage- 
ment recommendations for bioenergy crops are not 
always well defined. For example, nitrogen fertilizer 
recommendations for Miscanthus have been widely 
debated (Arundale et al., 2014a), with important impli- 
cations for the overall GHG footprint of production 
(Roth et al., 2015). To the extent that there are interac- 
tions between best management practices and site- 
level ecosystem properties (soil texture, land-use his- 
tory, etc.), assumptions about management used 
within an assessment study should ideally be imple- 
mented at the level of management decision-making, 
that is, the field-scale. Some landscape modeling stud- 
ies have started to endogenize management intensity 
questions by simulating productivity and environmen- 
tal impacts at different rates of N application (Gelfand 
et al., 2013) or different levels of tillage (Zhang et al., 
2010) for a given simulation run, for example, account- 
ing for realistic variations in best management prac- 
tices across the study landscape. 
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Study goals 


Our study used the DayCent biogeochemical process 
model to assess perennial grass productivity and associ- 
ated biogenic GHG emissions as a function of land qual- 
ity and management intensity. Implications for 
bioenergy landscape design were explored through a 
case study of switchgrass production around a newly 
constructed commercial-scale cellulosic biorefinery in an 
area with substantial heterogeneity in soils and land 
use. This investigation expands on previous work in 
two main ways: 

1. We conducted an extensive model parameterization 
and validation effort based on a data-diverse, multi- 
criteria approach, using a large parameterization 
dataset collected from the literature spanning wide 
gradients of climate, soil texture, and management 
intensity. 

2. We evaluated the impacts of management intensifica- 
tion at the full spatial resolution of the assessment, 
estimating optimal levels of nitrogen fertilizer appli- 
cation for each point on the landscape in order to 
either maximize yield or minimize biogenic GHG 
emissions. 


Our objective was to develop a rigorous, well-vali- 
dated spatially explicit biogeochemical modeling capa- 
bility that can serve as the basis for future integrated 
assessment and landscape optimization efforts. 


Materials and methods 


Case study introduction 


We performed a landscape assessment case study simulating 
the cultivation of switchgrass (Panicum virgatum) to supply bio- 
mass to the Abengoa cellulosic biorefinery located outside the 
town of Hugoton in southwestern Kansas (Peplow, 2014), 
which began operations in fall 2014. While the plant will ini- 
tially produce 25 million gallons of ethanol per year using corn 
stover as the primary lignocellulosic feedstock, switchgrass has 
been mentioned as an advanced cellulosic feedstock of interest 
and Biomass Crop Assistance Program Project Area 7 is spon- 
sored by the company and targets the establishment of 20 000 
acres of switchgrass production in the area (U.S. Department of 
Agriculture FSA, 2011). The case study region has long been at 
the center of agricultural sustainability and energy issues, hav- 
ing been deeply affected by the Dust Bowl in the 1930s (Kansas 
Historical Society, 2014) and being the site of the earliest 
hydraulic fracturing trials in the U.S. (Borowski, 2012). 

Today, the surrounding Stevens County is a highly diverse 
and productive agricultural area. In 2012, 21.4% of the county 
was dedicated to highly-productive irrigated corn cultivation 
(average yield of 12.1 Mg ha™', or 192 bushels acre~') and 11% 
to dryland wheat (1.1 Mg ha™!, or 18 bu acre~!), with smaller 
fractions devoted to other crops and pasture/rangeland, 


supporting an inventory of 45 500 head of cattle including 
calves (U.S. Department of Agriculture NASS, 2014). In this 
study we investigated the biogeochemical implications of con- 
verting non-irrigated cropland and rangeland in Stevens 
County and its six neighboring counties in southwestern Kan- 
sas and the Oklahoma panhandle to switchgrass cultivation 
(see Fig. 1), examining tradeoffs between productivity and 
associated biogenic GHG emissions as a function of underlying 
soil type and management intensity, specifically nitrogen fertil- 
izer application rate. Issues of land ownership, conservation 
easement status, and other land-use policy limitations are 
excluded here, but explored in a subsequent publication dedi- 
cated to bioenergy landscape optimization. 


The DayCent model 


Productivity and net fluxes of biogenic CO, and NO from soils 
under switchgrass cultivation were modeled with the DayCent 
biogeochemistry model (Parton ef al., 1998; Del Grosso et al., 
2011). DayCent is a semi-empirical process-based model that sim- 
ulates cycling of C, N, and water in natural and agroecosystems 
based on site-specific biophysical factors, land-use history, and 
management practices (e.g. tillage, fertilizer application, irrigation, 
etc.). The spatial and temporal scope of the model lies in between 
that of dedicated crop growth models (Miguez et al., 2012) and 
generalized earth climate system models (Anderson et al., 2013; 
Hallgren et al., 2013). DayCent has been used extensively to pre- 
dict yields and environmental impacts of switchgrass cultivation 
(Adler et al., 2007; Chamberlain et al., 2011; Davis et al., 2012; Lee 
et al., 2012) and is also used to predict agricultural soil GHG emis- 
sions for the annual Inventory of U.S. Greenhouse Gas Emissions 
and Sinks (U.S. Environmental Protection Agency, 2014a). 

DayCent computes soil temperature and moisture for differ- 
ent layers of the soil profile (resolved separately based on soil 
texture, bulk density, and pH) using daily climate data inputs. 
Crop growth (net primary productivity, or NPP) is simulated 
using species-specific parameters describing photosynthetic 
efficiency, tissue C/N ratio limits, above- and below-ground C 
partitioning, and phenology, many of which are determined 
empirically using model parameterization datasets as described 
previously. Daily biomass growth potential (NPPpot) is derived 
from top-of-atmosphere radiation (srad), corrected for atmo- 
spheric transmission losses and multiplied by a series of 0-1 
factors that represent deviations from ideal temperature (fiemp) 
or soil moisture (fmo) or limitations due to canopy immaturity 
or self-shading (fcanopy). This potential growth is then adjusted 
down if the available soil mineral nitrogen supply is limiting, 
as determined based on a maximum incremental biomass 
C : N ratio at the given level of plant maturity: 


NPP = max(NPPpot; NPPNiim) where: 


NPP ot = srad - transmission - ftemp * ft20 ` feanopy 


NPPNiim = Navail . (C/N) 


max 
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Fig. 1 Map of all switchgrass field trial sites included in the full model parameterization and validation dataset prior to ecotype/lati- 
tude filtering. Ring size indicates the number of experimental treatments (e.g. different ecotypes or nitrogen fertilizer application 
rates) conducted at that site, and color represents the types of data available. The 7-county case study area is highlighted in pink, with 


a star designating the biorefinery location. 


In addition to plant productivity DayCent also estimates soil 
carbon and nitrogen cycling including net changes in SOC 
levels and N2O emissions, the major constituents of the agricul- 
tural soil GHG balance. Carbon dynamics are simulated for soil 
surface pools and the top soil layer (0-20 cm by default), with 
organic matter represented by two litter pools (metabolic and 
structural) and three SOC pools- an ‘active’ pool representing 
microbial biomass and associated microbial products with a 
rapid turnover rates, and two others representing chemically/ 
physically stabilized carbon with decadal- (‘slow’ pool) and 
century-scale (‘passive’ pool) turnover times. Nitrogen mineral- 
ization and immobilization rates for each pool are controlled 
by the maximum and minimum permissible C : N ratios for 
each pool, soil temperature, soil moisture, and microbial effi- 
ciency as a function of soil texture (Parton et al., 1987). The soil 
nitrogen balance considers synthetic N fertilizer addition, man- 
ure and other organic N amendments, atmospheric deposition, 
volatilization, leaching, plant uptake, and N mineralization and 
immobilization associated with soil organic matter transforma- 
tions. The model simulates nitrification of ammonium (NH,*) 
to nitrate (NO3_) including NOx and NO by-products, as well 
as denitrification of nitrate to gaseous products (N20, NOx, 
Nə). Nitrification is simulated by multiplying available soil 
ammonium by a maximum potential nitrification rate adjusted 
based on soil temperature, water-filled pore space (WFPS), and 
pH limitations (Del Grosso et al., 2000; Parton et al., 2001). The 
overall rate of denitrification and the N2O/N;, ratio of its prod- 
ucts are modeled based on the availability of nitrate and 
organic matter substrates (as inferred from heterotrophic respi- 
ration rate), and local soil micropore redox state and gas 


diffusion rates as inferred from WFPS and heterotrophic respi- 
ration rate. 


Ecosystem model parameterization & validation 


We undertook an extensive parameterization and validation of 
the DayCent switchgrass growth model to improve overall sim- 
ulation accuracy and verify the model response across gradi- 
ents of land quality and management intensity. A large dataset 
of switchgrass field trials for the continental United States was 
assembled from the peer-reviewed literature. Studies were 
included in the dataset provided they (a) specified the underly- 
ing soil in sufficient detail that a corresponding Natural 
Resource Conservation Service Soil Survey Geographic data- 
base (SSURGO; Ernstrom & Lytle, 1993) map unit could be 
identified using the Web Soil Survey (U.S. Department of Agri- 
culture NRCS); (b) specified key management variables 
(switchgrass ecotype planted, N fertilizer application rate, har- 
vest date) in sufficient detail to define key DayCent simulation 
parameters; and (c) reported disaggregated annual yield and 
GHG results (studies reporting results averaged across multi- 
ple sites or N treatments were excluded). All modeled-vs.- 
measured performance comparisons were made on the basis of 
time-averaged results across the duration of the field trial, so 
studies that report either annual yield data or treatment- 
averaged data were included. If multiple cultivars of the same 
ecotype were included in a study, results were averaged across 
cultivars into a single representative value for the ecotype for 
simplicity. Additional details on the parameterization and 
validation dataset are provided in the Appendix S1. 
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A total of 24 appropriate studies were identified, covering 
573 annual biomass yield points across 147 unique combina- 
tions of site (soil, climate) and management treatments (N rate, 
harvest date, etc.). Initial exploratory analysis confirmed the 
need to exclude the first two seasons of yield data before 
switchgrass yields stabilize (Lesur et al., 2013; Arundale et al., 
2014b), and to filter out treatments for ecotypes grown outside 
their typical latitude range (up to 54°N for lowland varieties 
and down to 34°N for upland, as per Casler, 2012). The remain- 
ing dataset was then randomly split at the level of individual 
studies 70 : 30 into a parameterization dataset and an indepen- 
dent validation dataset (Table 1). Several of these field trials 
included more detailed information on biomass partitioning or 
nitrogen content (Frank ef al., 2004; Dohleman et al., 2012; 
Anderson-Teixeira et al., 2013), long-term changes in SOC 
based on either repeat measurements or paired sites (Ma et al., 
2001; Liebig et al., 2008; Follett et al., 2012; Anderson-Teixeira 
et al., 2013; Bonin & Lal, 2014), and/or time-resolved N,O 
emissions (Nikiema et al., 2011; Hong et al., 2012; Schmer et al., 
2012; Smith et al., 2013). There was insufficient data available to 
perform an independent validation of these model performance 
criteria, so all of these studies were included in the parameteri- 
zation dataset. While the raw SOC dataset was very noisy, 
eliminating individual data points that were not reported as 


statistically significant yielded a final reduced modeling dataset 
that behaved more predictably and was consistent with recent 
switchgrass SOC meta-analysis results (i.e. Qin et al., 2015a). 

The combined parameterization and validation dataset 
covers a wide range of latitudes and longitudes (Fig. 1) and 
temperature/precipitation regimes (Fig. S1), as well as a wide 
range of soil textures and Natural Resource Conservation Ser- 
vice land capability class (LCC) ratings (Helms, 1992; Fig. S2). 
Note that LCC ratings are reflective of a variety of land-use 
limitations, some of which are explicitly simulated in the Day- 
Cent model (e.g. dry climates, extreme textures, shallow soils) 
and some of which are not (e.g. erosion susceptibility, drainage 
class). The switchgrass crop parameterization was further 
informed with data from greenhouse or growth chamber exper- 
iments looking at productivity response across gradients of 
temperature (Balasko & Smith, 1971; Hsu et al., 1985; Reddy 
et al., 2008; Kandel et al., 2013; Wagle & Kakani, 2014) or mois- 
ture (Xu et al., 2006). 

The parameterization process started with a default switch- 
grass crop parameter set based on previous work (Adler et al., 
2007; Davis et al., 2012) and focused on refining parameters 
relating to productivity, temperature and moisture stress 
response, nitrogen management, shoot vs. root partitioning, and 
tissue death and turnover rates, with separate parameterizations 


Table 1 DayCent switchgrass crop parameterization and validation and data sources 


Parameter types Data sources 


Parameterization 
Aboveground biomass yield 


Ma et al. (2001), Fuentes & Taliaferro (2002), 


# treatments 


98 (67 after filtering) 


Vogel et al. (2002), Frank et al. (2004), Pearson (2004), 
Mulkey et al. (2006), Schmer et al. (2008, 2012), 
Nikiema et al. (2011), Dohleman et al. (2012), 

Follett et al. (2012), Hong et al. (2012), 

Kering et al. (2012), Anderson-Teixeira et al. (2013), 
Boyer et al. (2013), Bonin & Lal (2014) and 

Pedroso et al. (2014) 


Soil organic carbon (SOC) changes 


Ma et al. (2001), Liebig et al. (2008), Follett 


18 (8 after filtering) 


et al. (2012), Anderson-Teixeira et al. (2013) 
and Bonin & Lal (2014) 


Soil nitrous oxide (N20) emissions 


Nikiema et al. (2011), Hong et al. (2012), 11 


Schmer et al. (2012) and Smith et al. (2013) 


Seasonal aboveground and 
belowground biomass 
accumulation, and/or C/N ratio 
Phenology 


Frank et al. (2004), Dohleman et al. (2012) 3 
and Anderson-Teixeira et al. (2013) 


Sanderson (1992), Hopkins et al. (1995), NA 


Sanderson et al. (1997), Berdahl et al. (2005), 
Casler et al. (2007) and Wang et al. (2013) 


Productivity response to temperature 


Balasko & Smith (1971), Hsu et al. (1985), NA 


Reddy et al. (2008), Kandel et al. (2013) 
and Wagle & Kakani (2014) 


Productivity response to soil moisture 
Independent validation 


Yield Staley et al. (1991), Muir et al. (2001), 


Xu et al. (2006) 


NA 


49 (44 after filtering) 


Cassida et al. (2005), Adler et al. (2006), 
Fike et al. (2006), Arundale et al. (2014b) 
and Wilson et al. (2014) 
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for both upland and lowland ecotypes as appropriate. Initial 
exploratory analysis suggested that capturing differences in 
phenology between upland and lowland ecotypes was essential 
for accurate yield simulation, consistent with the current under- 
standing of maturation based on photoperiod being a strong 
determinant of yield differences between different cultivars 
grown at a given latitude (Casler et al., 2004). We set green-up 
dates uniformly for both ecotypes as function of latitude based 
on a variety of literature sources (Sanderson, 1992; Sanderson 
et al., 1997; Wang et al., 2013) as illustrated in Fig. S5. Peak bio- 
mass dates were predicted as a function of both latitude and 
ecotype as inferred from a variety of sources (Sanderson, 1992; 
Hopkins et al., 1995; Sanderson et al., 1997; Vogel et al., 2002; 
Frank et al., 2004; Berdahl et al., 2005; Casler et al., 2007; Ander- 
son-Teixeira et al., 2013; see Figs S6 and S7) and used to trigger 
plant senescence events within DayCent. 

After crop phenology was set other parameter adjustments 
were implemented manually, using the model automation rou- 
tine described in the next section to rapidly evaluate parameter 
changes against the 67-point parameterization dataset. The 
parameterization process focused on maximizing modeled-vs.- 
measured fidelity (based on visual inspection of the plotted 
data and calculation of modeled-vs-measured regression 
parameters, Pearson correlation coefficient, and root mean 
squared error; e.g. Smith et al., 1996) for upland and lowland 
ecotype yields, changes in SOC, and growing season N2O emis- 
sions, but also took into account time-resolved shoot : root and 
C : N ratio data where available for a multi-criteria evaluation 
of parameter set performance (e.g. Fig. $13). Once the parame- 
terization process was complete, independent validation of 
upland and lowland yield performance was conducted based 
on the data held in reserve (i.e. using the holdout validation 
method). 


Spatial data inputs, model initialization & automation 


A variety of spatially explicit data inputs are necessary to ini- 
tialize and run the DayCent model for a large-scale parameteri- 
zation or landscape analysis, including data on climate, soil 
type, and land-use history. Data sources used in this analysis 
are summarized in Table 2. Soil texture, rock fraction, and pH 


Table 2 Summary of spatial data inputs 


for different soil profile layers of the dominant soil component 
for each map unit were taken directly from the SSURGO data- 
base (Ernstrom & Lytle, 1993), and bulk density, field capacity, 
wilting point, and saturated hydraulic conductivity were com- 
puted using the Saxton equations (Saxton et al., 1986). Climate 
data on a 32 km grid was derived from the North American 
Regional Reanalysis database (NARR; Mesinger et al., 2006). 

Land-use history and current land management practices 
were compiled from a variety of sources. Current land use was 
determined from the National Land Cover Database 2006 
(NLCD; Wickham et al., 2013), re-sampled from the native 
30 m resolution to 240 m for ease of use, and re-classified into 
the simplified categories of annual agricultural lands (‘culti- 
vated crops’, ‘pasture/hay’), rangeland (‘dwarf shrub’, ‘shrub/ 
scrub’, ‘grassland/herbaceous’, ‘sedge/herbaceous’), and 
excluded (all other categories including forested and developed 
lands). Irrigated areas were identified using the MIrAD-US 
database (Pervez & Brown, 2010), and federally-owned lands 
were identified using the USGS Federal Lands of the United 
States data layer (U.S. Geological Survey, 2015) and excluded 
from further analysis (3.4% of the landscape, most part of 
Cimarron National Grassland). These five GIS layers were then 
intersected and small slivers were eliminated by merging all 
polygons smaller than 1 ha into the neighbor with which they 
shared the longest border. This yielded 39 320 polygons of a 
variety of sizes across the seven-county Hugoton case study 
area (Fig. S16), representing 3779 unique combinations of 
model inputs requiring individual simulation, which we refer 
to as DayCent modeling ‘strata’. 

For each strata, the DayCent model was pre-initialized using 
the same pre-settlement and historical land use assumptions as 
used in the EPA Inventory of U.S. Greenhouse Gas Emissions 
and Sinks (U.S. Environmental Protection Agency, 2014a) and 
described in detail by Ogle et al. (2010). Model initialization 
included an equilibration run of several thousand years dura- 
tion reflecting the natural state of the land prior to conversion 
to agriculture in order for all soil C and N pools to achieve 
steady-state values. Historical management between initial 
plow-out and the modern period was simulated with crop rota- 
tions and management practices compiled at regional scale 
from a variety of historical and modern sources (Ogle et al., 


Spatial Data 

database type Year Native resolution URL 

SSURGO Soils 2012 1 : 12 000-1 : 63 360 http:/ /www.nrcs.usda.gov/wps/portal/nrcs/detail/ 

soils/survey/geo/?cid=nres142p2_053627 

NARR Daily 1979-2010 32 km http://www.emc.ncep.noaa.gov/mmb/rreanl 
weather 

NLCD Land use 2006 30 m http://www.mrlc.gov/nlcd2006.php 

MIrAD-US Irrigation 2007 250 m http://earlywarning.usgs.gov/USirrigation/ 
extent 

Federal Lands Federalland 2005 640 acres/1 mi?/1 : 2 000 000 _ http://nationalmap.gov/small_scale/mld/fedlanp.html 

of the ownership 


United States 


MIrAD-US, MODIS Irrigated Agriculture Dataset for the United States; NARR, North American Regional Reanalysis; NLCD, National 


Land Cover Database; SSURGO, Soil Survey Geographic database. 
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2010). The future switchgrass simulations were then executed 
across part of a 29-node, 288-processor cluster computing sys- 
tem at the Colorado State University Natural Resource Ecology 
Laboratory. Parallel execution was implemented in Python 
(http://www.python.org/) using forking operations to take 
advantage of multiple cores within a given node. 


Landscape design analysis scenarios, results processing 
& sensitivity analysis 


For the landscape analysis case study we simulated conversion 
of all non-irrigated, non-federally owned polygons within the 
seven-county case study area to rain-fed lowland switchgrass 
cultivation. We conducted 30-year forward simulations to 
assess long-term productivity and trends in soil C and N 
cycling, recycling the full range of the NARR historic weather 
record to represent future weather conditions. In order to 
assess the response of crop productivity and GHG performance 
to management intensity, seven different rates of nitrogen fer- 
tilizer application were simulated for each strata (0- 
150 kgN ha~! in 25 kgN ha! increments). We assumed that 
switchgrass would be replanted every 10 years after field 
preparation consisting of chisel plow and field cultivator opera- 
tions, and that the crop would be neither fertilized nor har- 
vested the year of establishment in order to limit competition 
from weeds and ensure robust crop establishment, as per local 


Table 3 DayCent switchgrass crop parameter changes 


Type of parameter change 


Productivity potential and temperature response: 
Increase in productivity to reflect updated solar 
radiation model and compensate for increased 
belowground partitioning; differentiation in 
productivity and temperature response between 
upland and lowland ecotypes 

Growth response to moisture stress: Reduced 
sensitivity to soil moisture stress 

Belowground partitioning: Reduction in baseline 
BG partitioning rate. BG partitioning in response 
to moisture stress reduced, but response to 
nutrient stress increased 

Tissue N and lignin content: Root maximum 
allowable C : N ratio lowered slightly. Root lignin 
content reduced. Small amount of N fixation 
added to make growth under no-fertilizer 
conditions more realistic 

Tissue death rates: Death rate and fall rate for 
shoots increased. Root maturation rate increased, 
and turnover rate of both juvenile and mature 
roots increased 


N conservation: Increased translocation of 
nitrogen from shoots to roots during senescence 


extension recommendations. These assumptions are highly 
conservative as switchgrass is often established in this region 
without tillage, the need for periodic replanting is widely 
debated, and first-year harvest can be possible if the crop 
achieves sufficient first-year productivity. 

Switchgrass harvest yields, changes in SOC levels, and 
annual N20 emissions were then averaged over the 30-year sim- 
ulation period for each strata. Average annual NO emission 
values were converted into CO, equivalents using a 100-year 
global warming potential (GWP100) value of 298 (Forster et al., 
2007) then added to the CO) flux values associated with average 
annual net SOC changes for an estimate of total direct biogenic 
emissions. To determine the biogenic GHG intensity of produc- 
tion (Mg COzeq Mg”! biomass harvested) total emissions per 
hectare were divided by the associated simulated switchgrass 
yield. Continuous functions of biomass yield and biogenic GHG 
intensity vs. nitrogen fertilizer application rate were developed 
for each strata by applying a cubic regression to the 30-year 
averaged simulation results for the different N rate simulations, 
and the yield-maximizing and GHG balance-minimizing N 
rates were then interpolated. The simulated yield and GHG 
intensity associated with these optimized N rates was then re- 
associated with the appropriate landscape polygons, and aggre- 
gated across the full landscape in order to develop curves illus- 
trating total potential landscape productivity and biogenic GHG 
emissions balance when strata are managed for these different 


Parameter Original Lowland Upland 
name value value value 
PRDX(1) 2.75 4* 3.0% 
PPDF(1) 30 30 30 
PPDF(2) 45 44 44 
PPDF(3) 1 0.75 0.75 
PPDF(4) 2.5 2 2 
CWSCOEFKF(1,1) 0.38 0.35 0.35 
CWSCOEFK(1,2) 9 14 14 
CFRTCN(1) 0.5 0.7 0.7 
CFRTCN(2) 0.3 0.25 0.25 
CFRTCW(1) 0.6 0.4 0.4 
CFRTCW(2) 0.3 0.25 0.25 
PRBMX(1,1) 55 50 50 
FLIGNI(1,2) 0.26 0.06 0.06 
FLIGNI(1,3) 0.26 0.13 0.13 
SNFXMX(1) 0 0.005 0.005 
FSDETH(3) 0.05 0.075 0.075 
FALLRT 0.01 0.1 0.1 
CMXTURN 0.12 0.3 0.3 
RDRJ 0.4 0.72 0.72 
RDRM 0.2 0.54 0.54 
CRPRTF(1) 0.15 0.43 0.43 


*Change in the PRDX(1) values also reflects a recent change in the DayCent model to simulate atmospheric transmission losses of 


photosynthetically active radiation. 
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objectives. The sensitivity of these landscape results to key crop 
parameters, landscape characterization, and switchgrass culti- 
vation scenario assumptions was assessed as detailed in Table 4 
in order to determine the overall robustness of our conclusions. 
All results analysis routines were automated in Python through 
a combination of SQLite database operations (http: //sqlite.org/ 
) via the sqlite3 module (http://docs.python.org/2/library /sqli- 
te3.html), data manipulation in the native Python list data type, 
and figure generation using the matplotlib.pyplot module 
(http: //matplotlib.org/api/pyplot_api-html). 


Results 


Ecosystem model parameterization & validation 


A total of 79 different switchgrass parameter set itera- 
tions were ultimately developed and tested. Final 
upland and lowland ecotype parameter values that dif- 
fer from the default DayCent switchgrass crop parame- 
terization are detailed in Table 3. The most significant 
changes were: 

e Increased plant potential NPP rate — The most recent 
version of DayCent explicitly models solar radiation 
atmospheric transmission losses, and the revised 
potential NPP parameter PRDX(1) must be adjusted 
higher relative to previously-published versions to 
reflect the new growth calculation on a canopy pho- 
tosynthetically active radiation (PAR) basis rather 
than a top-of-atmosphere PAR basis. Further fine 
adjustments were made to PRDX(1) to optimize the 
observed yield difference between the different eco- 
types and to offset slightly increased belowground C 
partitioning. 

e Adjusted temperature and moisture stress response 
curves — Crop temperature (Fig. S8) and moisture 
stress (Fig. S9) response curves were set based on 
the greenhouse and growth chamber studies listed 
in Table 1. In the case of temperature response, 
fine adjustments to the edges of the curve where 
direct empirical data were lacking were imple- 
mented to improve overall modeled-vs.-measured 
yield performance across the full parameterization 
dataset. The same temperature and moisture 
response curves were used for both ecotypes. 
Comparison of measured and modeled yield 
ranges binned by site average growing degree day 
accumulation (Fig. S10) or annual precipitation 
(Fig. S11) verifies that the model accurately cap- 
tures increasing switchgrass productivity at war- 
mer, wetter sites. 


e Increased belowground partitioning and root turn- 
over — The default parameterization slightly under- 
estimated belowground biomass, significantly 
underestimated observed SOC increases, and over- 


predicted NO emissions. A small increase in 
belowground partitioning coupled with a large 
increase in root turnover rates resulted in more 
carbon being cycled into the soil and more mineral 
N being taken up by the plant, improving model 
performance on all three criteria. 


Model parameterization and validation results for 
yields and soil GHG fluxes are shown in Fig. 2. Suffi- 
cient data was available to perform holdout method 
independent validation of yield predictions for both the 
upland and lowland ecotypes (Fig. 2a). The out-of-sam- 
ple validation root mean square errors (RMSE) are 3.7 
and 4.1 Mg ha! for the upland and lowland ecotypes, 
respectively, with minimal bias (mean difference of 
—0.2 and +2.3 Mg ha“, respectively). When these data 
are binned by nitrogen fertilizer application rate 
(Fig. 2c) we see that the model is able to capture the 
general trend toward increased switchgrass productiv- 
ity with increasing management intensity. In contrast, 
yield response to land quality was more ambiguous. 
Neither measured nor modeled yields exhibited a clear 
relationship with soil texture across the full parameteri- 
zation and validation dataset (Fig. 2d), and the weak 
trend towards lower yields at sites with higher LCC rat- 
ing (more marginal) was not replicated in our model 
(Fig. S12). 

Field trial data on soil GHG balance were sparser, 
and the results presented in Fig. 2b reflect within-sam- 
ple model performance against the parameterization 
dataset rather than an independent validation. 
Observed annual changes in SOC under switchgrass 
were much larger than measured growing season N20 
emissions when compared in CO>-equivalent terms. The 
within-sample RMSE value for the combined SOC and 
NO dataset is 0.69 MgCO, eq ha! yr! across the 
wide range of sites, climates, and nitrogen application 
rates represented in the underlying studies (n = 19, 
Fig. 1). Bias calculations show that the model tends to 
err in the direction of overestimating N2O and underes- 
timating SOC accumulation, so the resulting predictions 
of switchgrass GHG balance are somewhat conserva- 
tive. Additional detail on SOC and NO performance is 
available in Figs 514 and S15. 


Landscape design case study 


Simulated lowland switchgrass yields as a function of 
nitrogen fertilizer application rate across the 3779 
unique DayCent strata in the Hugoton case study area 
are shown in Fig. 3 for both crop land and rangeland 
conversion, color-coded to the underlying soil texture 
of each strata. The maximum attainable yield under 
arbitrarily well-fertilized conditions shows significant 
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Fig. 2 Parameterization and validation results: (a) modeled vs. measured switchgrass biomass yield fits for lowland ecotype within 
sample (‘Low-param’, r = 0.10, root mean square error (RMSE) = 6.2 Mg ha’, mean difference (MD) = —2.1 Mg ha™') and out-of- 
sample [‘Low-valid’, r = 0.66, RMSE = 4.1, MD = +2.3] field trial results and upland ecotype within-sample (‘Up-param’, r = 0.38, 
RMSE = 4.0, MD = +2.1) and out-of-sample (‘Up-valid’, r = 0.26, RMSE = 3.7, MD = —0.2) results; (b) modeled vs. measured annual- 
ized changes in soil organic carbon (r = 0.08, RMSE = 1.05 MgCO2eq ha™! yr~', MD = +0.42 MgCO2eq ha™! yr~') and growing sea- 
son N20 emissions (r = 0.54, RMSE = 0.22, MD = +0.10) on a CO»-equivalent basis (combined r = 0.79, RMSE = 0.69, MD = +0.24); (c) 
measured and modeled switchgrass yield ranges for different levels of nitrogen application across all parameterization and validation 
data points; and (d) measured and modeled switchgrass yield ranges for different soil surface texture classes across all parameteriza- 


tion and validation data points. 


variation with soil texture, ranging from more than 
10 Mg ha“ in certain clay and sandy soils down to 
about 6 Mg ha! in the more moderate-textured silty 
soils. In this semi-arid climate of intermittent precipita- 
tion events soil moisture levels are often near wilting 
point, and simulated average yields reflect a tension 
between the greater total water holding capacity of 
finer-textured soils (an advantage during relatively 
wet years) vs. more effective infiltration and less 
surface soil evaporation in coarser soils (the so-called 
‘inverse texture effect’, beneficial during dry years; 


Noy-Meir, 1973; Epstein et al., 1997; Lane et al., 1998). 
Switchgrass yields on converted rangeland are gener- 
ally less sensitive to fertilizer application rates than 
those on converted cropland, in some cases showing 
no response to increasing N rate. This is due to higher 
background mineral nitrogen levels from transient soil 
organic matter turnover following conversion in these 
areas. However, for most strata, full switchgrass yield 
potential is realized at N rates from 60 to 
100 kgN ha™!, as indicated with circular markers in 
the upper panels of Fig. 3. 
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Fig. 3 Simulated yield and GHG intensity response to increasing nitrogen fertilizer rate for switchgrass production on former crop 
land and former rangeland across the 3779 distinct DayCent strata in the case study landscape. Interpolated yield maxima and GHG 
intensity minima are marked with solid markers for cropland conversion and open markers for rangeland conversion. The color of 
the lines and markers indicates the soil surface texture of the strata, with yellow = sand, red = clay, green = silt, and brown = loam. 


When examining landscape assessment results under 
different levels of management intensity it is impor- 
tant to note a discrepancy between the metrics of total 
area GHG emissions or mitigation (MgCOzeq per hec- 
tare per year) and the GHG intensity of biomass pro- 
duction (MgCO eq per Mg biomass grown). Biomass 
production GHG intensity results are shown as a func- 
tion of nitrogen application rate across all simulation 
strata in the lower panels of Fig. 3. For most strata at 
most nitrogen fertilizer rates, 30 year-average soil car- 
bon sequestration outweighs nitrous oxide emissions 
on a CO,-equivalent basis. Biomass has the lowest 
(most negative) direct biogenic GHG footprint when 
cultivated with no nitrogen fertilizer on previously 
cropped fine-textured soils due to their high potential 
for SOC accumulation. This sequestration value is 
increasingly offset at higher N rates due to marginal 
soil N,O emissions outpacing corresponding yield and 
SOC gains. In contrast, rangeland has higher initial 
SOC levels and thus less capacity for carbon seques- 
tration after conversion to bioenergy feedstock crop- 
ping, which results in an overall net GHG intensity 
much closer to zero. 

Total area biogenic GHG emissions are explored in 
relation to productivity in Fig. 4. Since switchgrass can 


be managed either to maximize yields or to optimize 
soil GHG balance (i.e. to maximize soil carbon seques- 
tration while minimizing NO emissions), we selected 
a random 10% of the landscape polygons, determined 
the nitrogen application rate that maximizes yield and 
the rate that maximizes GHG mitigation for each poly- 
gon, and aggregated the results for these different 
management strategies across our landscape sub-sam- 
ple. The difference between the two curves is repre- 
sentative of the degree to which system productivity 
and GHG performance can be modulated by adjusting 
management on the same limited land base. Convert- 
ing 10% of the landscape to switchgrass managed 
under either strategy results in >0.75 Mt biomass feed- 
stock production annually, enough material to supply 
approximately two facilities the size of the Abengoa 
biorefinery, with the associated net GHG impact of 
sequestering >0.06 Mt (60 000 metric tons) of CO». Fig- 
ure 5 illustrates that the optimal switchgrass manage- 
ment intensity for maximum GHG mitigation is 
related to the distribution of soil types and current 
land use across the case study area, factors that are 
also correlated with one another as conventional crop- 
ping tends to be concentrated on the moderate 
textured soils of this landscape. 
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Sensitivity analysis 


Our modeled landscape productivity and GHG balance 
are highly sensitive to certain key crop model parame- 
ters, particularly optimal growth temperature, potential 
NPP, and response to moisture stress (Table 4), high- 
lighting the importance of careful parameterization and 
validation. While the first two parameters were reason- 
ably well constrained as a result of the switchgrass 
parameterization process, empirical data on switchgrass 
moisture stress in the literature is rare, and uncertainty 
around those parameters in the model should be con- 
sidered large. Interestingly, the aggregate landscape 
simulation results show little sensitivity to the length of 
the switchgrass growing season in early spring and late 
fall, even though accurate characterization of the phe- 
nology proved essential for determining crop parameter 
sets with adequate out-of-sample validation perfor- 
mance. Accurate characterization of the landscape in 
terms of climate and land-use history is important to 
the integrity of simulation results as well, with GHG 
mitigation particularly sensitive to past land-use his- 
tory. Landscape results were generally less sensitive to 
crop management scenario assumptions around tillage 
intensity and fertilizer timing. 


Discussion 


Challenges of model development and data-diverse 
parameterization 


This study was grounded in an extensive model param- 
eterization and validation effort using a data-diverse 
multi-criteria approach enabled by automation of all 
simulation runs and results analyses. The diversity of 
studies included in our switchgrass field trial dataset 
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Fig. 4 Cumulative total soil GHG mitigation vs. cumulative 
switchgrass biomass production for a random 10% of the case 
study landscape under different management goals: managing 
each land parcel to maximize switchgrass yields (dashed line, 
highest-productivity sites aggregated first), or managing each 
land parcel to maximize ecosystem GHG mitigation (solid line, 
strongest mitigation sites aggregated first). 


was intended to provide a highly general test of model 
performance independent of a single study, environ- 
ment, or management practice, and the simultaneous 
consideration of different types of data (Table 1) was 
designed to ensure that accuracy of one model perfor- 
mance criteria was never improved at the expense of 
others. This approach proved useful for sorting out 
interrelated model responses, for example, the strong 
effects of belowground partitioning and root turnover 
rate parameters on harvested biomass yield, SOC 
changes, and N20 rates that must be balanced in a sys- 
tematic manner during model parameterization. 

Our resulting model explains approximately half of 
the observed variability in yield and GHG observations 
in our dataset, and realistically captures climate effects 
and responses to management intensity. While soil type 
has been observed to have an effect on bioenergy grass 
productivity in semi-arid climates outside the U.S. 
(Di Virgilio et al., 2007; o Di Nasso et al., 2015; Roncucci 
et al., 2015), we did not observe a strong texture signal 
in our domestic switchgrass productivity dataset, con- 
sistent with a similar previous large-scale model param- 
eterization effort (Wullschleger et al., 2010). Future 
modeling work on soil-climate interactions in the U.S. 
would greatly benefit from additional field trails like 
Wilson et al. (2014) that include paired trails across mul- 
tiple soil types on a landscape, in order to tease apart 
soil-precipitation interactions without the other larger 
sources of variability present in composite datasets. 

While we believe the approach presented here repre- 
sents an improvement over more limited model param- 
eterization efforts, it is still possible to over- 
parameterize to a dataset of this size and achieve good 
fits via unrealistic mechanisms that do not translate well 
out-of-sample. In our experience, reliable yield perfor- 
mance was only achieved once crop phenology was 
adequately captured. This parameterization dependency 
is somewhat challenging, as there are a limited number 
of studies in the literature reporting detailed phenology 
in the form required for this type of generalized model- 
ing effort (i.e. dates for green-up and peak biomass or 
senescence as a function of crop cultivar and site lati- 
tude). Additionally, since maturation dates can vary sig- 
nificantly even for cultivars within a single ecotype 
(Frank et al., 2004), focusing on only two ecotype group- 
ings introduces additional errors. Accurate representa- 
tion of phenology can be superseded to a certain extent 
by narrowing crop temperature response curves to trun- 
cate early and late-season productivity, but our initial 
parameterization efforts in this direction performed 
very poorly on validation. 

After several dozen parameter set iterations we were 
approaching the limit of what could be accomplished 
with setting parameter values manually, as model 
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Fig. 5 Map illustrating (a) soil texture distribution, (b) current land use, and (c) the switchgrass nitrogen fertilizer application rates 
associated with maximum area GHG mitigation (as opposed to biomass GHG intensity explored in Fig. 3), for all non-federal non-irri- 
gated cropland and rangeland across the Hugoton case study area. Aerial imagery is included as a background layer for scale, and 
the boundary of the 7-county case study area is shown in white with the biorefinery location marked with a yellow star. 


responses became more subtle and antagonistic across 
the different performance criteria. This, coupled with 
the relatively large number of crop parameters in the 
DayCent model, suggests that future efforts could bene- 
fit greatly from systematic parameter optimization tech- 
niques such as inverse modeling and/or x-fold cross 
validation. Such techniques provide a more systematic 
and transparent approach, facilitate maximum extrac- 
tion of information from a field trial dataset, and help to 
identify and avoid over-parameterization issues. While 
inverse modeling techniques have been demonstrated to 
improve model performance against small datasets 
(Necpalova et al., 2015), they are not commonly applied 
to crop growth and biogeochemistry models used in 
bioenergy sustainability analysis at the scale of a full 
parameterization and validation effort such as the one 
presented in this study. 

More broadly, model performance is dependent on 
underlying model structure in addition to parameteriza- 
tion, and DayCent model development efforts to 
improve its performance in perennial grass systems are 
ongoing. For example, accurate yield estimates depend 
on capturing crop phenology and canopy dynamics, 
and efforts are underway to represent early-season leaf 
canopy closure and late-season senescence dynamics 
more realistically. Similarly, while past model develop- 
ment efforts have largely focused on the upper 0-20 cm 
soil layer most relevant for annual crops, representation 
of soil carbon changes under perennial grasses may be 
improved through more explicit simulation of deep 


rooting, dissolved organic carbon movement, and other 
aspects of soil organic matter dynamics deeper in the 
soil profile (Campbell, 2015). 


Landscape design case study interpretation and climate 
accounting issues 


The Hugoton case study site was selected because it fea- 
tured the most heterogeneity in soils and current land 
use among the first three commercial-scale cellulosic 
biorefinery sites in the U.S. (the other two being located 
in the Iowa corn belt; Peplow, 2014). However, the dry 
climate at this site proved challenging from a modeling 
perspective, with landscape simulation results showing 
high sensitivity to interactions between crop moisture 
stress parameterization and soil texture, and with N20 
response to increasing N rate lower than what might be 
expected in a wetter climate. Overall, our landscape 
simulation results are similar to others in the literature 
that find the GHG balance of perennial grasses domi- 
nated by soil carbon sequestration immediately post- 
establishment, with only modest N20 emissions (e.g. 
Gelfand et al., 2013). In field trials where baseline N O 
emission rates are reported and an Intergovernmental 
Panel on Climate Change (IPCC) Tier 1 estimate can be 
made, observed switchgrass NO emissions rates are 
often near or below the lower bound of the IPCC NO 
emissions factor range (0.3— 3% fertilizer N emitted as 
N2O-N; Fig. S15). This suggests that switchgrass is an 
efficient nitrogen cycler, and that the critique of biofuel 
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Table 4 Sensitivity analysis 


Parameter or assumption Change 


Baseline average landscape performance 
Yield = 7.24 Mg ha™! yr! 
GHG mitigation = 0.76 Mg COveq ha“! yr7! 
DayCent crop parameterization 
Optimal temperature 
Productivity potential 
Moisture sensitivity 


Increase PPDF(1) by 3 degrees 
Decrease PRDX(1) value by 10% 
Switch to the default 


curve associated with annual crops 


Root death rates Decrease RDRJ and 
RDRM by 10% 
Decrease growing season 


length by 15 days 


Phenology 


Site characterization 
Land-use history Assume a uniform cropped 

history rather than a mix 

of cropped and grazed 

Assume the NARR weather 

for Hugoton for the entire 


7-county study area 


Climate 


Agronomic management 
Fertilizer date Assume fertilizer application 
1 month earlier 


Tillage intensity 


GHG, greenhouse gas; NARR, North American Regional Reanalysis. 


GHG mitigation benefits being outweighed by direct 
and indirect N?O emissions (Crutzen et al., 2008) is 
likely overstated for second-generation perennial grass 
feedstocks. Our predicted SOC sequestration rates after 
cropland conversion are very similar to median meta- 
analysis results in Qin et al. (2015a), while the finding of 
a small net positive sequestration with the conversion of 
grassland is consistent with, but on the optimistic edge 
of, their observed range. 

Full quantification of the confidence interval around 
these landscape assessment results is highly challenging 
due to multiple levels of uncertainty in problem specifi- 
cation, landscape input data, model structure, and 
model parameterization (Walker et al., 2003). While pre- 
vious work has shown that uncertainty in DayCent-esti- 
mated changes in SOC is dominated by model 
structural uncertainty (Ogle et al., 2010), repeating this 
empirical uncertainty assessment approach here is 
impractical as only one of the 145 treatments in our 
parameterization and validation dataset includes data 
for all three of the factors that determine biomass GHG 
intensity (yield, SOC change, and N20; Fig. 1). Here we 
rely on sensitivity analysis to identify priorities for 
future model improvement efforts. When discussing 
model sensitivity and uncertainty issues it is also impor- 
tant to be cognizant of issues in the underlying field 


Switch to no-till crop establishment 


Change in total 
landscape productivity 


Change in total landscape 
GHG mitigation 


—21% —42% 
—IAN —18% 
—1.1% +8.6% 
+0.11% +0.60% 
+0.49% —0.41% 
+5.0% +29% 
—4.2 —9.6% 
+0.11% +1.5% 
+0.13% +0.97% 


trial datasets. Field measurements of N2O are very chal- 
lenging due to emissions variability on extremely fine 
spatial (Li et al., 2013) and temporal (Jorgensen et al., 
2012; van der Weerden et al., 2013) scales, and studies 
based on sampling at weekly or every-other-week fre- 
quency (as was the case for all N20 studies in our data- 
set) are vulnerable to systematic biases of up to 20% 
and 60%, respectively (Parkin, 2008; Smemo et al., 2011). 
Additionally, small-scale field trials do not necessarily 
reflect imperfections in agronomic management (e.g. 
uneven fertilizer application) and real-world harvest 
losses, possibly introducing a systematic over-estima- 
tion of switchgrass productivity (Searle & Malins, 2014). 

Finally, there are several issues around climate impact 
accounting relevant to this landscape study. While soil 
carbon sequestration will eventually attenuate even 
though N20 emissions will persist for as long as N fer- 
tilizer is being applied (Sheehan et al., 2003; Adler et al., 
2007), our assessment for this semi-arid system suggests 
it will take 60-80 years for annual sequestration and 
N20 rates to reach parity (Fig. S17). A more dynamic 
climate impact accounting approach that takes transient 
forcing benefits into account (e.g. Holtsmark, 2015) 
would tend to further weight near-term SOC benefits 
against future N»O emissions, though no standardized 
accounting approach is yet widely accepted (U.S. Envi- 
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ronmental Protection Agency, 2014b). The current study 
is also somewhat limited in focusing on the climate 
impacts of biomass feedstock production from a purely 
biogeochemical perspective, ignoring potential biophysi- 
cal impacts such as changes in surface albedo or evapo- 
transpiration and water dynamics that are significant in 
some bioenergy production scenarios (Munoz etal., 
2010; Georgescu et al., 2011; Cherubini et al., 2012; 
Caiazzo et al., 2014), and ignoring any potential broader 
impacts on other ecosystem services (Chamberlain & 
Miller, 2012). Additionally, future changes to atmo- 
spheric CO, concentrations and climate are not consid- 
ered here, though they may have large repercussions 
for landscape design (Bryan et al., 2010). The ability of 
current ecosystem models to accurately extrapolate to 
such future conditions is an active area of investigation 
(De Kauwe et al., 2013). Additional sensitivity analysis 
on these points would improve our understanding of 
bioenergy landscape performance potential, and high- 
light priorities for future research efforts. 
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